Reference: Forecasting: Principles and Practice
Exploring & Visualizing Time Series
Scotty Technologies Inc. (“Scotty”) adalah sebuah perusahaan start-up teknologi yang didirikan pada tahun 2017 di Istanbul, Turkey. Salah satu layanan utama dari Scotty adalah motorcycle ride-sharing atau yang akrab kita ketahui Ojek Online. Menurut informasi yang dilansir dari markets.businessinsider.com, Scotty berencana untuk menjadi super-app pertama di Turkey dengan tekonologi dan model bisnis yang distruptif, serta berencana untuk melanjutkan pertumbuhan yang menjanjikan dan mencapai profitabilitas dalam waktu dekat. Berbicara mengenai pertumbuhan, tentu saja tidak luput dari demand. Oleh karena itu, projek ini ditujukan untuk membuat model analisa dan memprediksi demand transaksi per-jam pada layanan motorcycle ride-sharing di Scotty. Dalam proses pembuatan model ini, kita menggunakan data transaksi Scotty dari 2017-10-01 sampai 2017-12-02.
Problem yang diselesaikan yaitu membuat model untuk memprediksi demand transaksi per-jam pada layanan motorcycle ride-sharing di Scotty. Namun yang hendak diprediksi dalam rentang waktu berapa lama? Mari kita cek dulu supaya dapat membantu proses pemodelan yang akan kita buat.
#> Observations: 504
#> Variables: 3
#> $ src_sub_area <chr> "sxk97", "sxk97", "sxk97", "sxk97", "sxk97", "sxk97", ...
#> $ datetime <dttm> 2017-12-03 00:00:00, 2017-12-03 01:00:00, 2017-12-03 ...
#> $ demand <lgl> NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA...
Struktut data diatas belum sesuai. Mari sesuaikan struktur datanya dan meilhat rentang waktu yang hendak diprediksi dahulu:
scotty_test <- scotty_test %>%
mutate(
# Transaksi Per-Jam
datetime = floor_date(datetime, unit="hour"),
demand = as.integer(demand)
)
summary(scotty_test)#> src_sub_area datetime demand
#> Length:504 Min. :2017-12-03 00:00:00 Min. : NA
#> Class :character 1st Qu.:2017-12-04 17:45:00 1st Qu.: NA
#> Mode :character Median :2017-12-06 11:30:00 Median : NA
#> Mean :2017-12-06 11:30:00 Mean :NaN
#> 3rd Qu.:2017-12-08 05:15:00 3rd Qu.: NA
#> Max. :2017-12-09 23:00:00 Max. : NA
#> NA's :504
Data diatas merupakan summary dari data yang hendak diprediksi. Berdasarkan informasi tersebut, dapat disimpullkan bawah masalah yang perlu kita selesaikan yaitu melakukan prediksi demand per-jam di sub area sxk97, sxk9e dan sxk9 dalam rentang waktu dari 2017-12-03 00:00:00 sampai 2017-12-09 23:00:00 atau 24 jam * 1 minggu. Oke, mari kita mulai proses nya.
Berikut ini adalah data yang kita punya untuk melakukan analisa dan prediksi.
5 Top Line Data
5 Bottom Line Data
| Variable | Description |
|---|---|
| id | Transaction id |
| trip_id | Trip id |
| driver_id | Driver id |
| rider_id | Rider id |
| start_time | Transaction Time |
| src_lat | Request source latitude |
| src_lon | Request source longitude |
| src_area | Request source area |
| src_sub_area | Request source sub-area |
| dest_lat | Requested destination latitude |
| dest_lon | Requested destination longitude |
| dest_area | Requested destination area |
| dest_sub_area | Requested destination sub-area |
| distance | Trip distance (in KM) |
| status | Trip status (all status considered as a demand) |
| confirmed_time_sec | Time different from request to confirmed (in seconds) |
Missing Value Jika dilihat terdapat missing value pada variabel trip_id dan driver_id, namun tidak menjadi masalah karena pada case ini kita butuh data waktu, sub area dan demand.
#> id trip_id driver_id rider_id
#> 0 4676 4676 0
#> start_time src_lat src_lon src_area
#> 0 0 0 0
#> src_sub_area dest_lat dest_lon dest_area
#> 0 0 0 0
#> dest_sub_area distance status confirmed_time_sec
#> 0 0 0 0
Duplicate Value Dataset ini tidak memiliki data transaksi yang duplikat, sehingga bisa kita lanjutkan ke tahap pre-processing.
data.frame(jumlah.observasi.awal=length(scotty_input$id),
jumlah.observasi.unik=length(unique(scotty_input$id)))#> Observations: 90,113
#> Variables: 16
#> $ id <chr> "59d005e1ffcfa261708ce9cd", "59d0066affcfa261708...
#> $ trip_id <chr> "59d005e9cb564761a8fe5d3e", NA, "59d006c131e39c6...
#> $ driver_id <chr> "59a892c5568be44b2734f276", NA, "599dc0dfa5b4fd5...
#> $ rider_id <chr> "59ad2d6efba75a581666b506", "59cd704bcf482f6ce2f...
#> $ start_time <dttm> 2017-10-01 00:00:17, 2017-10-01 00:02:34, 2017-...
#> $ src_lat <dbl> 41.07047, 41.07487, 41.04995, 41.05287, 41.06760...
#> $ src_lon <dbl> 29.01945, 28.99528, 29.03107, 28.99522, 28.98827...
#> $ src_area <chr> "sxk9", "sxk9", "sxk9", "sxk9", "sxk9", "sxk9", ...
#> $ src_sub_area <chr> "sxk9s", "sxk9e", "sxk9s", "sxk9e", "sxk9e", "sx...
#> $ dest_lat <dbl> 41.11716, 41.08351, 41.04495, 41.08140, 41.02125...
#> $ dest_lon <dbl> 29.03650, 29.00228, 28.98192, 28.98197, 29.11316...
#> $ dest_area <chr> "sxk9", "sxk9", "sxk9", "sxk9", "sxk9", "sxk9", ...
#> $ dest_sub_area <chr> "sxk9u", "sxk9e", "sxk9e", "sxk9e", "sxk9q", "sx...
#> $ distance <dbl> 5.379250, 1.126098, 4.169492, 3.358296, 11.69357...
#> $ status <chr> "confirmed", "nodrivers", "confirmed", "confirme...
#> $ confirmed_time_sec <dbl> 8, 0, 32, 65, 0, 27, 32, 30, 0, 24, 9, 0, 118, 1...
Dataset yang dimiliki terdiri dari 90.133 observasi dan 16 variabel. Jika dilihat dari struktur datanya, beberapa variabel memiliki tipe data yang belum sesuai. Namun, karena projek ini bertujuan untuk melakukan time-series forecasting terhadap demand per-jam pada layanan motorcycle ride-sharing di Scotty, maka data yang kita perlukan adalah waktu transaksi per-jam dan demand pada setiap src_sub_area. Maka dari Berikut prosesnya:
Dalam melakukan analisa dan prediksi time series ada beberapa hal yang harus terpenuhi antara lain:
1. Data tidak boleh ada yang missing
2. Data harus terurut berdasarkan periode waktunya
3. Tidak boleh ada waktu atau periode yang terlewat/bolong
Terkait case ini, Jam opearasional Scotty yaitu 24jam/hari, sehingga kita harus memastikan seluruh data jam sudah lengkap secara sekuensial dari rentang tanggal minimum hingga tanggal maximum. Jika terdapat deret waktu yang kosong, maka dapat dilakukan imputasi data dengan asumsi bahwa pada waktu tersebut memang tidak terdapat transaksi, sehingga jumlah demand = 0.
scotty_input <- scotty_input %>%
group_by(src_sub_area) %>%
# PAD tanggal sekuensial
padr::pad(start_val = min(scotty_input$datetime), end_val = max(scotty_input$datetime)) %>%
ungroup() %>%
distinct() %>%
mutate(
# replace NA Value menjadi 0
demand = ifelse(is.na(demand),0,demand)
) %>%
arrange(datetime)Berikut data hari yang tidak memiliki transaksi:
Data diatas menunjukan total waktu yang tidak terdapat transaksi dari sub Area sxk97, sxk9e dan sxk9s yaitu 301 Jam.
#> src_sub_area datetime demand
#> Length:4536 Min. :2017-10-01 00:00:00 Min. : 0.00
#> Class :character 1st Qu.:2017-10-16 17:45:00 1st Qu.: 5.00
#> Mode :character Median :2017-11-01 11:30:00 Median : 15.00
#> Mean :2017-11-01 11:30:00 Mean : 19.87
#> 3rd Qu.:2017-11-17 05:15:00 3rd Qu.: 27.00
#> Max. :2017-12-02 23:00:00 Max. :217.00
Informasi di atas merupakan summary keseluruhan data yang kita miliki. Dapat diketahui bahwa dataset ini merupakan data transaksi dari 3 sub area yaitu sxk97, sxk9e dan sxk9s yang terjadi dari 2017-10-01 00:00 sampai 2017-12-02 23:00 dengan jumlah minimum 0 demand dan jumlah maximum 217 demand. Demand 0 diasumsikan bahwa pada waktu tertentu memang tidak terdapat transaksi.
Berikut ini visualisasi demand berdasarkan sub-area dari dari 2017-10-01 00:00 sampai 2017-12-02 23:00. Titik/Poin berwarna merah, hijau dan biru menandakan data deret waktu tersebut kosong atau tidak terdapat transaksi pada waktu tersebut.
scotty_demand <- scotty_input
mean_demand <- scotty_demand %>%
group_by(src_sub_area) %>%
summarise(mean = mean(demand))
scotty_demand <- sqldf("SELECT SD.*, MD.mean
FROM scotty_demand SD
INNER JOIN mean_demand MD ON SD.src_sub_area = MD.src_sub_area")
scotty_demand <- scotty_demand %>%
mutate(
weekdays = weekdays(datetime),
popup = glue(
"Category: {src_sub_area}
Transaction Time: {datetime}
Day: {weekdays}
Demand: {demand}"
))
plot_scotty_demand <- ggplot(scotty_demand,aes(x=datetime, y=demand),show.legend = FALSE)+
#geom_point()+
geom_line(show.legend = FALSE,size=0.3)+
geom_point(aes(text=popup),size=0.05)+
geom_point(data = scotty_demand %>% filter(demand==0),
aes(x=datetime, y=demand, color=src_sub_area, text=popup),size=0.5, show.legend = FALSE)+
labs(
x = NULL,
y = NULL,
title = "Demand by Sub-Area"
) +
scale_x_datetime(date_breaks = "7 day", labels = date_format("%b %d"))+
facet_wrap(~ src_sub_area, scale = "free", ncol = 1)+
theme_minimal()+
theme(
legend.position = "none",
plot.title = element_text(hjust = 0.5),
panel.spacing = unit(1, "lines"))
ggplotly(p = plot_scotty_demand, tooltip="text") %>%
layout(showlegend=FALSE)Pada grafik diatas dapat dilihat bahwa jumlah demand dari ketiga sub-area relatif sama dan juga terdapat pola seasonal meskipun belum terlihat jelas detailnya. Selain itu, juga dapat dilihat demand pada Sabtu,4 November 2017 pukul 01.00 dan Sabtu, 25 November 2017 pukul 01.00 jauh lebih tinggi dari hari lainnya. Sayangnya data yang kita miliki belum dapat menjawab hal tersebut, namun mari kita fokus ke demand. Berikut adalah demand per-harinya:
# Daily Demand
polar_day <- scotty_input %>%
mutate(day=weekdays(datetime)) %>%
group_by(src_sub_area,day) %>%
summarise(demand_per_day=sum(demand)) %>%
ungroup() %>%
mutate(
day = ordered(day, levels=c("Monday","Tuesday","Wednesday","Thursday","Friday","Saturday","Sunday"))
)
serror <- function(x) sqrt(var(x)/length(x))
ggplot(polar_day %>% filter(src_sub_area=="sxk97"), aes(x=day, y=demand_per_day, fill = day)) +
geom_bar(width = 1, stat = "identity", color = "white", show.legend = FALSE) +
geom_errorbar(aes(ymin = demand_per_day - serror(demand_per_day),
ymax = demand_per_day + serror(demand_per_day),
color = day),
width = .2) +
scale_y_continuous(breaks = 0:nlevels(day)) +
labs(
title = "Daily Demand in Sub-Area sxk97"
)+
theme_minimal() +
theme(axis.title = element_blank(),
legend.position = "none",
plot.title = element_text(hjust = 0.5, size=12, face="bold"),
axis.text.y = element_blank(),
axis.text.x=element_text(face="bold"))+
coord_polar() -> day_polar1
ggplot(polar_day %>% filter(src_sub_area=="sxk9e"), aes(x=day, y=demand_per_day, fill = day)) +
geom_bar(width = 1, stat = "identity", color = "white", show.legend = FALSE) +
geom_errorbar(aes(ymin = demand_per_day - serror(demand_per_day),
ymax = demand_per_day + serror(demand_per_day),
color = day),
width = .2) +
scale_y_continuous(breaks = 0:nlevels(day)) +
labs(
title = "Daily Demand in Sub-Area sxk9e"
)+
theme_minimal() +
theme(axis.title = element_blank(),
legend.position = "none",
plot.title = element_text(hjust = 0.5, size=12, face="bold"),
axis.text.y = element_blank(),
axis.text.x=element_text(face="bold"))+
coord_polar() ->day_polar2
ggplot(polar_day %>% filter(src_sub_area=="sxk9s"), aes(x=day, y=demand_per_day, fill = day)) +
geom_bar(width = 1, stat = "identity", color = "white", show.legend = FALSE) +
geom_errorbar(aes(ymin = demand_per_day - serror(demand_per_day),
ymax = demand_per_day + serror(demand_per_day),
color = day),
width = .2) +
scale_y_continuous(breaks = 0:nlevels(day)) +
labs(
title = "Daily Demand in Sub-Area sxk9s"
)+
theme_minimal() +
theme(axis.title = element_blank(),
legend.position = "none",
plot.title = element_text(hjust = 0.5,size=12, face="bold"),
axis.text.y = element_blank(),
axis.text.x=element_text(face="bold"))+
coord_polar() ->day_polar3
grid.arrange(day_polar1,day_polar2,day_polar3, ncol = 3)Grafik diatas memberikan informasi total demmand per-hari dari masing-masing Sub-Area. Secara keseluruhan, ketiga Sub-Area menunjukan pola demand per-hari yang sama, dimana total demand paling besar yaitu pada hari Jumat dan total demand paling kecil yaitu pada hari Minggu. Mari kita lihat data demand per-jamnya berikut:
# Hourly Demand
polar_hourly <- scotty_input %>%
mutate(hour=hour(datetime)) %>%
group_by(src_sub_area,hour) %>%
summarise(demand_per_hour=sum(demand)) %>%
ungroup() %>%
mutate(
hour = ordered(hour,levels=c("0","1","2","3","4","5","6","7","8","9","10","11","12","13","14","15","16","17",
"18","19","20","21","22","23"))
)
ggplot(polar_hourly %>% filter(src_sub_area=="sxk97"), aes(x=hour, y=demand_per_hour, fill = hour)) +
geom_bar(width = 1, stat = "identity", color = "white", show.legend = FALSE) +
geom_errorbar(aes(ymin = demand_per_hour - serror(demand_per_hour),
ymax = demand_per_hour + serror(demand_per_hour),
color = hour),
width = .2) +
scale_y_continuous(breaks = 0:nlevels(hour)) +
labs(
title = "Hourly Demand in Sub-Area sxk97",
subtitle = "24 Hour Format"
)+
theme_minimal() +
theme(axis.title = element_blank(),
legend.position = "none",
plot.title = element_text(hjust = 0.5,size=12, face="bold"),
plot.subtitle = element_text(hjust = 0.5,size=11),
axis.text.y = element_blank(),
axis.text.x=element_text(size=11, face="bold"))+
coord_polar() -> hour_polar1
ggplot(polar_hourly %>% filter(src_sub_area=="sxk9e"), aes(x=hour, y=demand_per_hour, fill = hour)) +
geom_bar(width = 1, stat = "identity", color = "white", show.legend = FALSE) +
geom_errorbar(aes(ymin = demand_per_hour - serror(demand_per_hour),
ymax = demand_per_hour + serror(demand_per_hour),
color = hour),
width = .2) +
scale_y_continuous(breaks = 0:nlevels(hour)) +
labs(
title = "Hourly Demand in Sub-Area sxk9e",
subtitle = "24 Hour Format"
)+
theme_minimal() +
theme(axis.title = element_blank(),
legend.position = "none",
plot.title = element_text(hjust = 0.5,size=12, face="bold"),
plot.subtitle = element_text(hjust = 0.5,size=11),
axis.text.y = element_blank(),
axis.text.x=element_text(size=11, face="bold"))+
coord_polar() -> hour_polar2
ggplot(polar_hourly %>% filter(src_sub_area=="sxk9s"), aes(x=hour, y=demand_per_hour, fill = hour)) +
geom_bar(width = 1, stat = "identity", color = "white", show.legend = FALSE) +
geom_errorbar(aes(ymin = demand_per_hour - serror(demand_per_hour),
ymax = demand_per_hour + serror(demand_per_hour),
color = hour),
width = .2) +
scale_y_continuous(breaks = 0:nlevels(hour)) +
labs(
title = "Hourly Demand in Sub-Area sxk9s",
subtitle = "24 Hour Format"
)+
theme_minimal() +
theme(axis.title = element_blank(),
legend.position = "none",
plot.title = element_text(hjust = 0.5,size=12, face="bold"),
plot.subtitle = element_text(hjust = 0.5,size=11),
axis.text.y = element_blank(),
axis.text.x=element_text(size=11, face="bold"))+
coord_polar() -> hour_polar3
grid.arrange(hour_polar1,hour_polar2,hour_polar3, ncol = 3)Grafik diatas memberikan informasi total demmand per-jam dari masing-masing Sub-Area. Secara keseluruhan, ketiga Sub-Area menunjukan pola demand per-jam yang sama, dimana total demand relatif tinggi dari pukul 17:00 hingga pukul 19:00 dan total demmand relatif rendah dari pukul 01:00 hingga pukul 0:00. Namun, yang berbeda adalah total demand sub-area sxk9e juga cukup tinggi pada pukul 08:00, hal ini berbeda dengan sub-area lainnya.
Sebelum melakukan forecasting kita perlu mengetahui informasi tren dan seasonal yang ada pada data. Salah satu cara untuk mengetahuinya yaitu melakukan Decompose. Decompose merupakan tahapan dalam analisis time series yang digunakan untuk menguraikan beberapa komponen yang berupa:
Jika pada hasil decompose, trend masih membentuk sebuah pola maka dapat dicurigai masih ada seasonality yang belum ditangkap oleh frekuensi data. Seharusnya trend cenderung naik atau cendurung turun atau membentuk garis yang smooth.
Sesuai dengan masalah yang hendak diselesaikan yaitu melakukan prediksi demand per-jam pada masing-masing sub-area selama 1 minggu, maka mari kita coba decompose secara harian dan mingguan terhadap salah satu sub-area.
sxk9s <- scotty_input %>% filter(src_sub_area == "sxk9s") %>% select(-src_sub_area)
sxk9s_ts_weekly <- sxk9s %>% .$demand %>% ts(frequency = 24*7)
plot_sxk9s_ts_weekly <- sxk9s_ts_weekly %>%
tail(24*7*4) %>%
decompose() %>%
autoplot()+
theme_minimal()+
theme(
title = element_text(size=11),
axis.title=element_text(size=9, face="bold"),
axis.text.x=element_text(size=9),
axis.text.y=element_text(size=9),
legend.position = "none",
panel.spacing = unit(4, "pt")
)
ggplotly(plot_sxk9s_ts_weekly) %>%
layout(title = list(text = paste0('Decomposition on Weekly in Sub-Area sxk9s',
'<br>',
'<sup>',
'Deompose single seasonal time-series object using frequency = 24*7 (Weekly)',
'</sup>')))Hasil decompose objek single sesonal time-series diatas menggunakan frekuensi Mingguan pada data 1 bulan terkahir. Beberapa informasi yang bisa kita dapat yaitu:
Hasil decompose diatas mengindikasikan bahwa data time-series yang kita observasi memiliki pola Multi-Seasonal. Berdasarkan pola seasonalnya menunjukan bahwa data ini memiliki seasonal harian dan mingguan, sehingga mari kita mebuat data multiple seasonal time series dan melakukan decompose kembali.
sxk9s_msts_weekly <- sxk9s %>% .$demand %>%
msts(seasonal.periods = c(24, 24*7))
plot_sxk9s_msts_weekly <- sxk9s_msts_weekly %>%
tail(24*7*4) %>%
mstl() %>%
autoplot()+
labs(
title = "Decomposition on Daily and Weekly in Sub-Area sxk9s"
)+
theme_minimal()+
theme(
title = element_text(size=12),
axis.title=element_text(size=9, face="bold"),
axis.text.x=element_text(size=9),
axis.text.y=element_text(size=9),
legend.position = "none",
panel.spacing = unit(4, "pt")
)
ggplotly(plot_sxk9s_msts_weekly) %>%
layout(title = list(text = paste0('Decomposition on Daily & Weekly in Sub-Area sxk9s',
'<br>',
'<sup>',
'Deompose multiple seasonal time-series object using seasonal.periods = c(24, 24*7) (Daily & Weekly)',
'</sup>')))Pada hasil decompose objek multiple seasonal time-series diatas, dapat dilihat pola tren terlihat lebih smooth dan tidak membentuk pola. Hal ini mengindikasikan bahwa kemungkinan besar seluruhan seasonal sudah ditangkap. Untuk lebih jelasnya, berikut ini visualisasi seasonal perjam-nya dalam harian dan mingguan pada data observasi kita:
data_msts <- msts(scotty_input %>% filter(src_sub_area=="sxk9s") %>%
select(demand), seasonal.periods = c(24, 24*7))
sxk9s_msts_weekly_dc <- mstl(data_msts)
plot_seasonal_msts <- as.data.frame(sxk9s_msts_weekly_dc) %>%
mutate(
datetime = sxk9s$datetime
) %>%
mutate(
day = wday(datetime, label = TRUE, abbr = FALSE),
hour = as.factor(hour(datetime))
) %>%
group_by(day, hour) %>%
summarise(
seasonal = sum(Seasonal24 + Seasonal168)
) %>%
ggplot(mapping = aes(x = hour, y = seasonal)) +
geom_col(aes(fill = day)) +
scale_fill_viridis_d(option = "plasma") +
labs(
title = "Seasonal Analysis based Daily and Weekly in Sub-Area sxk9s",
fill = "Day of Week"
) +
theme_minimal() +
theme(
title = element_text(size=11),
axis.title=element_text(size=9, face="bold"),
axis.text.x=element_text(size=9),
axis.text.y=element_text(size=9),
legend.position = "right"
)
ggplotly(plot_seasonal_msts)